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Artificial fragmentation of the matter density field causes the formation of spurious groups of 
particles in N-body simulations of non-standard Dark Matter (DM) models which are characterized 
by a small scale cut-off in the linear matter power spectrum. These spurious halos alter the prediction 
of the mass function in a range of masses where differences among DM models are most relevant 
to observational tests. Using a suite of high resolution simulations we show that the contamination 
of artificial groups of particles significantly affect the statistics of halo spin, shape and virial state 
parameters. We find that spurious halos have systematically larger spin values, are highly elliptical 
or prolate and significantly deviate from virial equilibrium. These characteristics allow us to detect 
the presence of spurious halos even in non-standard DM models for which the low-mass end of 
the mass function remains well behaved. We show that selecting halos near the virial equilibrium 
provides a simple and effective method to remove the bulk of spurious halos from numerical halo 
catalogs and consistently recover the halo mass function at low masses. 


I. INTRODUCTION 

A generic feature of dark matter models alternative to 
the standard Cold Dark Matter (CDM) hypothesis is the 
presence of a small scale cut-off in the linear matter power 
spectrum. In the Warm Dark Matter (WDM) scenario 
this arises from the damping of matter density fluctu¬ 
ations on scales smaller than the free-streaming length 
of WDM particles (see e.g. PQ ;2j). Similarly, in mod¬ 
els of Ultra-Light Axion (ULA) dark matter, the effec¬ 
tive Jean’s length associated to the oscillating axion field 
naturally introduces a damping scale in the matter power 
spectrum (see e.g. [3]). Scenarios of Late Forming Dark 
Matter (LFDM) also predict a suppression of power at 
small scales that depends on the epoch of DM formation 
[3]. Because of this feature, non-standard DM models 
and the standard CDM scenario differ in the late-time 
clustering properties at small scales, where the dynam¬ 
ics of the gravitational collapse is highly non-linear. In 
this regime quantitative predictions have been possible 
thanks to the use of large volume high-resolution N-body 
simulations (see e.g. [SHE]). 

The accuracy of N-body results relies upon the abil¬ 
ity to control numerical systematic effects. These can be 
controlled through convergence analysis tests that eval¬ 
uate the dependence on the volume and mass resolution 
of the simulations. However, in the case of models with a 
sharp cut-off in the initial power spectrum, N-body sim¬ 
ulations have shown the presence of unphysical group of 
particles which result from the artificial fragmentation of 
the matter density field at small scales [THUTH] , These 
artifacts have been shown to form independently of how 
initial conditions are generated [TS]. This is because arti¬ 
ficial fragmentation is a discretization phenomenon which 
arises from N-body sampling Poisson noise at wave num¬ 
bers larger than the cut-off in matter power spectrum. 

Artificial halos contaminate the low-mass end of the 
halo mass function in the interval range where differ¬ 
ence among the different DM scenarios are most relevant. 
Hence, their removal is essential to predict the correct 


abundance of dwarf-galaxy halos as well as the dynami¬ 
cal properties of halo sub-structures that in recent years 
have become a probe of the nature of DM through ob¬ 
servations of the local Universe [TM251 . 

Numerical studies have shown that increasing the mass 
resolution of the simulations alleviates artificial fragmen¬ 
tation mm- I n fact, increasing the number of N-body 
particles (N p ) reduces the intra-particle distance d (i.e. 
increases the Nyquist frequency of the simulations) as 
well as the amplitude of the Poisson noise (oc 1 /N p ) thus 
leading to less fragmentation near the cut-off scale. How¬ 
ever, it would require the fluid limit to completely re¬ 
move this effect (for an attempt in this direction see e.g. 
ES -25). Hence, in the presence of artificial halos we 
ought to rely on empirical methods to infer artifact-free 
model predictions. 

Wang & White [18] have suggested to apply a mass- 
cut to numerical halo catalogs and retain only halos with 
mass larger than Mn m = 10.1pdfc“ e 2 ak , where p is the 
mean cosmic matter density and fc pea k is the wavenumber 
of the peak in the dimensionless power spectrum A 2 (fc) of 
WDM models and the numerical coefficient is estimated 
from simulations. 

Schneider, Smith and Reed [23] outlined an alterna¬ 
tive approach based on the fact that artificial halos con¬ 
tribute to an upturn of the low-mass end of the mass 
function with a characteristic power law behavior. Using 
this information, the contribution of artificial halos can 
be extrapolated at larger masses and subtracted from the 
measured mass function. 

A more sophisticated approach was adopted by Lovell 
et al. m- Their method combines multiple criteria: 
first, Lagrangian patches in the initial conditions are 
identified for each halo in the catalogs, those associated 
to patches whose shape is flatter than a given thresh¬ 
old are removed since genuine proto-halos are spheroidal; 
second, halos below a given mass-cut are removed; fi¬ 
nally, a match between the residual halos in simulations 
at different resolution is applied and only halos present at 
both resolutions are retained. Schneider m performed 
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a study of the properties of halos in non-standard DM 
scenarios using a similar approach. 

Here, we develop a complementary method that relies 
on the structural properties of halos as a way to differ¬ 
entiate genuine halos from artificial ones. Using a suite 
of high resolution N-body simulations of non-standard 
DM models we show that artificial halos are character¬ 
ized by extreme values of the spin and shape parameters 
and significantly depart from virial equilibrium. Hence, 
discarding unrelaxed halos from numerical halo catalogs 
can account for most of the effects induced by artificial 
group of particles and allows to predict the correct halo 
mass function over a wider range of mass than using sim¬ 
ple mass cuts. 

The paper is organized as follows: in Section [H] we 
present the models and describe the N-body simulation 
characteristics. In Section m we show the results on 
the halo mass function and the analysis of the structural 
properties of halos. In Section [IV] we propose a halo 
selection criterion based on the virial state of halos and 
present the results on the halo mass function and the dis¬ 
tribution of halo spins. Finally, we discuss and conclude 
in Section m 


II. N-BODY SIMULATIONS 

A. Models and Simulation Characteristics 

We run a suite of cosmological simulations using RAM¬ 
SES [28] , an adaptive mesh refinement code with a tree- 
based data structure that allows recursive grid refine¬ 
ment on a cell-by-cell basis. Particles are evolved using 
a particle-mesh (PM) solver, while the Poisson equation 
is solved using a multigrid method [25]. 

We consider as reference model a standard cold dark 
matter scenario with cosmological constant (ACDM) 
with zero curvature and the following set of parame¬ 
ters: n m = 0.3, h = 0.7, erg = 0.8, n s = 0.96 and 
f \ = 0.046. For the non-standard DM cosmologies we 
consider two WDM models with thermal relic particle 
mass of rn w( im = 1.465 keV (WDM-a) and 0.696 keV 
(WDM-b), and a LFDM model with phase transition 
redshift Zt = 1.5 x 10 6 (for details on the LFDM model, 
see HU). The cosmological parameters for these non¬ 
standard DM models are set to our reference ACDM val¬ 
ues. We note that out of these three non-standard DM 
models, only the LFDM satisfies the constraints imposed 
on the matter power spectrum inferred from the Lyman- 
a measurements m- WDM models with m w dm < 2 keV 
are disfavored at 4 a confidence level, however we include 
them in this paper to study artificial fragmentation of 
the matter density field. The linear matter power spec¬ 
tra of the non-standard DM models are characterized by 



FIG. 1: Linear matter power spectra at 2 = 0 for ACDM (red 
solid line), LFDM (blue dotted line), WDM-a (green dashed 
line) and WDM-b (brown dash-dotted line). 


a cut-off scale and damping slope 1 . These are shown for 
z = 0 in Fig. |T] We can see that the WDM-b model has 
a larger scale cut-off compared to WDM-a and LFDM, 
while the damping slope of LFDM is less steep than that 
of the WDM models. 

We use the linear power spectra to generate Gaussian 
initial conditions using the Zel’dovich approximation as 
implemented in MPGRAFIC [33] ■ To facilitate compar¬ 
ison we have used the same phase of the initial condi¬ 
tions for all simulated models. We set the starting red- 
shift of the simulation by imposing the constraint rela¬ 
tion o-(A™ arse ) = 0.02, such that at the initial redshift 
all models have the same value of the standard devia¬ 
tion of the initial density field a smoothed on the scale of 
the coarse grid A™ arse . This guarantees that the initial 
redshift of the simulations is sufficiently high to suppress 


1 We compute the initial linear power spectra using CAMB m • 
In the case of the WDM models we evaluate the linear transfer 
function as in m- 

T(k) = [l + {ak) 2u ] ~ " , (1) 

with spectrum cut-off scale 

a ^ m (^) 0 W (^ Y 3 (^) 115 (^) 0 - 29 ( 2 ) 

Vo Aj Vo.65/ Vm wdm ; Vswdm/ 

in units of h —1 Mpc. We assume damping slope v = 1 and 
number of degrees of freedom <7 w dm = 1.5. The power spectrum 
of WDM models is then obtained from that of ACDM as 

Pw dmW = T(/c) 2 Pacdm(&)- (3) 
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Model 

L (h -1 Mpc) 

Ap 

m p (h 1 M 0 ) 

ACDM 

27.5 

512 3 

1.29 x 10 7 

ACDM 

» 

1024 3 

1.61 x 10® 

LFDM 

55 

512 3 

1.29 x 10 7 

LFDM 

55 

1024 3 

1.61 x 10® 

WDM-a 

55 

512 3 

1.29 x 10 7 

WDM-a 

55 

1024 3 

1.61 x 10® 

WDM-b 

55 

512 3 

1.29 x 10 7 

WDM-b 

55 

1024 3 

1.61 x 10® 

WDM-b 

64 

512 3 

1.63 x 10 s 


velocity of the halo center of mass). The definition of 
the spin parameter given by Eq. Q normalizes the halo 
angular momentum to its maximal value. 

We also quantify the shape of halos by computing the 
symmetric mass distribution tensor defined as (see e.g. 

m) 





i=i 


r u,c)( r 0,i - rp ,c), 


( 6 ) 


TABLE I: N-body simulation characteristics. L is the simu¬ 
lation box length, N v is the number of N-body particles and 
m p the mass resolution. 


spurious effects due to transients |5ij . 

In Table [I] we list the characteristics of the N-body 
simulation runs. For all models we run (27.5 hu 1 Mpc ) 3 
volume simulations with 512 3 and 1024 3 particles corre¬ 
sponding to a mass resolution of m p = 1.29 x 10 7 h _1 Mg 
and uip = 1.61 x 10 6 h _ 1 M 0 respectively. For WDM-b 
we have also run a lower resolution simulation in larger 
volume: (64 h _1 Mpc ) 3 volume with 512 3 particles cor¬ 
responding to a mass resolution m p = 1.63 x 10 8 h^Mg. 
This simulation suite allows us to perform a detailed con¬ 
vergence study on artificial halos. 


B. Halo Finder and Halo Properties 


We detect halos using the code pFOF 35: based on a 
parallelized version of the friend-of-friend algorithm TT(ij . 
This detects halos as group of particles characterized by 
an intra-particle distance smaller than a given linking 
length parameter 6 , which we set to 0.2. During the first 
iteration, the algorithm groups all particles in the simula¬ 
tion box which are within distance b of an initial particle. 
This search is then carried out for the rest of the particles 
of the group until no new neighbors are found. This final 
set of particles is then tagged as a halo and removed from 
the list. The algorithm then iterates over the next un¬ 
tagged particles repeating the above procedure to detect 
a new halo. 

For each halo in the simulations we estimate the spin 
parameter defined as (551 : 


V2MVR' 


( 4 ) 


where M is the mass of the halo, V = \JGM/R is the 
orbital velocity at the virial radius R (G being Newton’s 
constant) and J is the amplitude of the halo angular 
momentum: 


N h 

J = m p '^2(f i - r c ) x (vi - v c ), (5) 

i=1 


where the sum runs over all IV^-halo particles with posi¬ 
tion ft and velocity Vi (r c and v c being the position and 


where a,/3 = 1,2,3 denote the three components of the 
position vectors. The eigenvalues of M a p, a 2 > b 2 > c 2 , 
define a triaxial ellipsoid with axis lengths a > b > c. 
These can be combined to define halo shape parame¬ 
ters. Here, we focus on the sphericity s = c/a , ellip- 
ticity e = ^(a — c)/(a + b + c) and prolateness p = 
\{a — 2b + c)/(a+ b + c). Genuine halos, being close to 
spheroidal, have large values of sphericity and low ellip- 
ticity and prolateness. On the other hand, spurious halos 
are expected to be abnormally elliptical (large value of 
e) with disk-like (oblate, p < 0 ) or needle-like (prolate, 
p > 0) shapes. 

We estimate the dynamical state of halos by computing 
the virial state parameter rj = 2K/\E\ £ [0, oo], where 
K is the total kinetic energy and E is the gravitational 
potential energy. In the case of virialized halos 77 « 1. It 
is worth reminding that 77 only provides an approximate 
estimation of a halo’s virial state. In fact, while there is 
no ambiguity in the determination of the kinetic energy 
of the halo (i.e. the sum of the kinetic energy of each halo 
particle), the potential energy E is a non-local quantity, 
since it also depends on particles that do not belong to 
the halo, but are in the surrounding density field. This 
may introduce a systematic source of uncertainty in the 
evaluation of the virial state. Despite such a limitation, 
p remains a useful proxy especially to identify halos with 
large deviation from the virial condition (i.e. 77 1 ). 


III. ANALYSIS OF SPURIOUS HALOS 

A visual inspection of the simulations shows that there 
are progressively less structures in LFDM, WDM-a and 
WDM-b compared to the ACDM. This can be seen in 
Fig-H where we plot a 2D projection of the dark matter 
density field at 2 = 0 in a box of ~ 7 h _1 Mpc side from 
the (27.5 h _1 Mpc ) 3 volume simulations with 512 3 parti¬ 
cles. This trend is consistent with expectations due to the 
fact that the damping of the linear matter power spec¬ 
trum of the non-standard DM models considered here 
occurs at an increasingly larger scale (small wavenum¬ 
bers) as we go from LFDM to WDM-b (see Fig. [lj. For 
these models we can also notice an increasing granularity 
of the density held along filamentary structures with the 
appearance of regularly spaced group of particles that are 
especially distinguishable in the WDM-b case. These are 
spurious numerical halos. 
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FIG. 2: 2D projection of the dark matter density field at z = 0 in a box of ~ 7 h x Mpc side from the simulations with (27.5 
h -1 Mpc) 3 volume and 512 3 particles for ACDM, LFDM, WDM-a and WDM-b (from left to right). 


A. Halo Mass Function 

The presence of a cut-off in the linear matter power 
spectrum introduces a characteristic mass scale in the 
halo mass function. As shown in [37] for WDM models, 
this scale can be parameterized in terms of the ‘half¬ 
mode’ mass Mhm, namely the mass associated to the 
length scale where the linear transfer function drops to 
half of its value. In the case of the WDM models consid¬ 
ered here we have Mh m = 7-41 x 10 9 h -1 M 0 for WDM-a 
and Mh m = 9.65 x 10 10 h -1 M 0 for WDM-b, while for 
the LFDM model we estimate Mh m = 3 x 10 9 h” 1 Mg. 

In Fig.[3]we plot the halo mass function at z = 0 for ha¬ 
los with at least 100 particles from the simulation suite. 
Let us first consider the results of the (27.5 h” 1 Mpc) 3 
volume runs with 512 3 particles. We can see that the 
number density of halos in the non-standard DM models 
is suppressed with respect to the ACDM case at around 
the half-mode mass as expected. The larger the half¬ 
mode mass the larger the suppression of the mass func¬ 
tion relative to ACDM. However, in the case of the WDM 
models we can see that the low mass end of the mass func¬ 
tion rises instead of falling. The effect is more dramatic 
for WDM-b model which has the largest half-mode mass. 
As shown in [13] this upturn is an artifact due to the con¬ 
tribution of artificial halos. Above the half-mode mass 
all models converge to the ACDM mass function within 
finite volume errors. Because of the absence of an up¬ 
turn in the LFDM case one might think that the model 
is exempt from spurious halo contamination. However, 
the mass function from the higher resolution run shown 
in Fig. [3] reveals this not to be the case since an upturn 
is present at much smaller masses. 

Let us compare the mass function from the (27.5 h -1 
Mpc) 3 volume simulations with 512 3 and 1024 3 particles. 
From Fig. [ 3 ] we see that the slope of the upturn at the 
low-mass end depends on the mass resolution of the sim¬ 
ulations, in agreement with the findings of [23]. This 
dependence is a clear indication of the fact that the ma¬ 
jority of halos in this mass range are numerical artifacts. 
We find that the mass functions of the lower resolution 
runs converge to those of higher resolution at ~ 5% for 


M> 4 x 10 9 h _1 M 0 - roughly corresponding to halos 
with more than 300 particles. For the WDM-b simu¬ 
lation of (64 h _1 Mpc) 3 volume and 512 3 particles we 
find a similar level of convergence for M> 5 x 10 ln h _1 
Mg, again roughly corresponding to a minimum of 300 
particles per halo. Thus, from now on we only consider 
halos with at least 300 particles. For simulations of (27.5 
h _1 Mpc) 3 volume and 512 3 particles, a 300-particle halo 
has a mass of 4 x 10 9 h _1 Mg and well below the half¬ 
mode mass Mhm of the WDM models (green and blue 
vertical dotted lines in Fig. [ 3 ]). The WDM mass func¬ 
tions, instead of falling, continue to rise even below their 
respective Mh m values. As such, even a conservative 300- 
particle cut does not fully solve the problem of spurious 
halos. Moreover, from the analysis of the LFDM mass 
function with 512 3 particles we can deduce that for a 
given mass resolution simulation the absence of a well 
defined upturn at the low-mass end of the mass function 
does not imply that a model characterized by a suppres¬ 
sion of power at small scales is exempt from artificial halo 
contamination. 

The point that we want to stress is that artifacts caused 
by spurious halos cannot be addressed by solely assuming 
conservative mass cuts in the halo catalogs. Moreover, 
extrapolating information from the upturn of the mass 
function might not be possible since depending on the 
cosmological model and the simulation characteristics an 
upturn may be absent near the half-mode mass. A case 
in point is the LFDM mass function (see Fig. [ 3 ]) from the 
(27.5 h -1 Mpc) 3 volume simulations with 512 3 particles, 
where at the lowest mass end there is no well defined 
upturn. As we will show hereafter, removing the contri¬ 
bution of artificial halos requires their characterization in 
terms of physical properties that may distinguish them 
from genuine halos. 


B. Halo Spin and Shape 

In Fig. [4] we plot the normalized probability den¬ 
sity function of the halo spin parameter at z = 0 
from the boxlength L = 27.5h _1 Mpc simulations 




5 



M [h- 1 MJ 

FIG. 3: Halo mass function at z = 0 for ACDM (red circles), 
LFDM (blue circles), WDM-a (green circles) and WDM-b 
(brown circles/stars) from the simulations listed in Table [I] 
Errorbars are given by Poisson errors. Vertical dotted lines 
indicate the value of the half-mode mass for LFDM (blue), 
WDM-a (green) and WDM-b (brown) models. 

with 512 3 particles. The different lines correspond to 
equally spaced logarithmic mass bins in the range 4 < 
M[10 9 IT 1 M 0 ] < 8 for ACDM, LFDM and WDM-a, and 
4 < M[10 9 h -1 M e ] < 100 for WDM-b 2 

Previous studies have shown that the distribution of 
halo spins in ACDM model depends only weakly on the 
halo mass and is approximately described by a log-normal 
distribution mm- As we can see from the top left panel 
of Fig. [4] our results are consistent with these findings. 
Notice that the tail of the ACDM distribution appears 
to be heavier in the lowest mass bin (red solid line). As 
shown in m this is mostly due to unrelaxed halos. We 
will confirm this conclusion in the analysis presented in 
Section IIII Cl 

In the top right panel of Fig. [4] we plot the distribution 
of spins for the LFDM model. Compared to the ACDM 
case, we can see a stronger dependence on halo mass. 
The distribution appears to be well described by a log¬ 
normal only for the highest mass bin (blue solid line), 
while becoming increasingly heavy tailed at lower masses. 

Despite a larger statistical noise, such a trend can also 
be seen in the WDM models shown in the lower pan¬ 
els of Fig. [4] We may notice large departures from log- 
normality. As in the LFDM case the distribution is in¬ 
creasingly heavy tailed, in particular the lower mass bins 


2 The upper limit of the mass range roughly corresponds to the 

largest half-mode mass of the pair of DM models considered. 



FIG. 4: Probability density function of the spin parameter at 
z = 0 for ACDM (top left panel), LFDM (top right panel), 
WDM-a (bottom left panel) and WDM-b (bottom right 
panel). The different lines correspond to 5 equally spaced 
logarithmic mass bins in the range 4 < M[lO 9 h - 1 M 0 ] < 
8 for the LFDM, ACDM and WDM-a models and 4 < 
M[10 9 h -1 M 0 ] < 100 for the WDM-b model. The red (blue) 
solid line corresponds to the lowest (highest) mass bin. 



- 0.1 0 0.1 0 . 2 - 0.1 0 0.1 0 . 20.1 0 0.1 0 . 20.1 0 0.1 0.2 


FIG. 5: Probability density function of halo sphericity, ellip- 
ticity and prolateness (panels top to bottom) at z = 0 for 
ACDM, LFDM, WDM-a and WDM-b (panels left to right) 
respectively. We consider the same mass bins as in Fig. |4] 

are characterized by a bimodal distribution of the spin 
parameter. Since spurious halos dominate the halo abun¬ 
dance in this mass range (see Fig. [3]), we can deduce that 
the bimodality is a numerical artifact similar to the up¬ 
turn in the halo mass function. Furthermore, as the low¬ 
est mass bin in the LFDM case in Fig. [4] shows a heavy 
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FIG. 6: Normalized probability density function in the A ' — r) plane at 2 = 0 for halos with mass M< 10 10 h 1 M@ (left panels) 
and 10 11 < M[1i -1 Mq] < 10 13 (right panels) in the case of ACDM, LFDM, WDM-a and WDM-b models. 


tail compared to a log-normal distribution, this indicates 
that spurious halos also affect the LFDM halo mass func¬ 
tion (see Fig. [ 3 ]) as confirmed by the presence of an upturn 
in the higher resolution run. This implies that other non¬ 
standard DM scenarios characterized by a suppression of 
power at small scales similar to LFDM, such as a class of 
Self-Interacting Dark Matter (SIDM) models discussed 
in H32> may be contaminated by spurious halos as well. 

In Fig. [5] we plot the normalized probability density 
function of halo shape parameters at 2 = 0 as a function 
of halo mass. We can see a mass trend similar to that 
of the spin parameter shown in Fig. [4] In the ACDM 
case the majority of halos are slightly elliptical, more¬ 
over the distribution of the shape parameters is nearly 
independent of the halo mass. This is not the case for 
the low mass bins of the non-standard DM models. In 
particular, in the WDM case we see that only the higher 
mass bins have a distribution that is consistent with that 
of the ACDM. At lower masses, halos have systemati¬ 
cally higher values of ellipticity and prolateness. This 
suggests that besides high spin values, artificial halos are 
highly elliptical and prolate compared to genuine ones. 
Though not directly comparable, these results are con¬ 
sistent with the analysis by Lovell et al. m which found 
that in WDM simulations spurious subhalos are associ¬ 
ated to proto-halos characterized by small values of the 
sphericity, while genuine subhalos at higher masses are 
more spheroidal and closer to that of ACDM prediction. 
This trend is similar to that of the sphericity distribution 
shown in the top panel of Fig. [5j Similar results were 
found in (27: which studied the distribution of shape pa¬ 
rameters in non-standard DM scenarios with power sup¬ 
pression at small scales. 


C. Virial State vs Spin Statistics 

Further insight on the structural properties of artifi¬ 
cial halos can be gained by considering the probability 
distribution of the virial state parameter 77. Since spu¬ 
rious halos are the result of artificial fragmentation, it 
is reasonable to expect that such group of particles are 
not virialized. Therefore, in the case of artificial halos 
the parameter 77 deviates from unity and correlates with 
large values of halo spin. 

In Fig. [ 6 ] we show the density plot of the normal¬ 
ized probability density function in the A' — 77 plane 
for the simulations with 512 3 particles and boxlength 
L = 27.5h _1 Mpc 2 = 0. The set of four panels on the 
left shows the distribution of halos with mass M< 10 10 
h _ 1 MQ, while the set on the right shows halos with mass 
in the range 10 11 < M[h _ 1 M Q ] < 10 13 . In the former 
case we can clearly see a strong correlation between large 
deviations from virial equilibrium (i.e. 77 = 1 ) and high 
spin values. This manifests through a heavy tail that is 
absent in the higher mass bin case. From the density 
plot we see that the fraction of halos populating this ex¬ 
tended tail is lowest for ACDM, while it increases for the 
non-CDM models in proportion to their half-mode mass. 
Moreover, the tail systematically extends towards larger 
values of A' and 77. We also notice that in the case of 
the WDM models the tail of the distribution becomes bi- 
modal. Let us remark that while the ACDM halos with 
M< 10 10 1i _ 1 Mq shown Fig. [b] have a distribution of 
spin and shape parameters which is nearly independent 
of their mass (see Figs. [4|5]), this is not the case for the 
non-standard DM models. This suggests that while the 
extended tail of the ACDM probability density function 
in the A 7 — 77 plane is composed of genuine unrelaxed halos, 
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that of the non-CDM models is dominated by spurious 
halos. 


IV. SPURIOUS HALOS SELECTION 

In the previous Section we have shown that the struc¬ 
tural properties of spurious halos in non-standard DM 
simulations differ from that of genuine halos 3 . In partic¬ 
ular, we have found a strong correlation between devia¬ 
tions from the virial condition and large values of halo 
spin at low masses. Spurious halos with such features 
contribute to the heavy tale of the spin distribution in 
the LFDM model and the bimodality in WDM models. 
These results suggest a simple way of removing artifi¬ 
cial halos from N-body halo catalogs since a cut based 
on the virial state of halos, as parametrized by 77 , can 
filter out spurious halos in the tail of the halo distribu¬ 
tion. In the CDM case this approach leaves only halos 
that are approximately well relaxed (see e.g. EDI). while 
it effectively removes the bulk of spurious halos in non- 
CDM models. As we will show here, this halo selection 
allows us to recover the halo mass function consistent 
with expectations at low masses and gives convergent re¬ 
sults when applied to N-body runs with different mass 
resolution. 



M [h-‘ M g 


3 Here, we have not considered the density profile of halos, since 
spurious groups of particles are irregularly shaped and unrelaxed. 
As shown in m the density profile of such perturbed systems 
is badly fit by the universal Navarro-Frenk-White profile, conse¬ 
quently the concentration parameter is completely uninformative 
about the internal structure of the halos. 


FIG. 7: Halo mass function at z = 0 for ACDM (red circles), 
LFDM (blue circles), WDM-a (green circles) and WDM-b 
(brown circles/stars) from the simulations listed in Table [T] 
before and after removing unrelaxed halos (grey points). All 
halos have at least 300 particles. Vertical lines correspond to 
the mass cut Mi; m from JT8] for WDM-a (green) and WDM-b 
(brown) models in the case of the simulations with boxlength 
L = 27.5h -1 Mpc and 512 3 particles. 

In Fig. [7] we plot the halo mass function at z = 0 
for halos with at least 300 particles for the simulations 
listed in Table [T] The grey points show the mass func¬ 
tions for halos with the virial state parameter in the 
range 0 < 77 < 1.5. In particular, let us focus on the 
mass function of LFDM (blue points), WDM-a (green 
points) and WDM-b (brown points) models before and 
after the virial state parameter cut from the simulations 
of boxlength L = 27.5 h _1 Mpc with 512 3 (empty circles) 
and 1024 3 (solid circles) particles. Notice that discard¬ 
ing unrelaxed halos (77 > 1.5) has removed the upturn at 
low masses. Now, the mass functions drop as expected 
below the half-mode mass of each model. Moreover, we 
can see that for the LFDM and WDM-a models, in the 
mass range where the halo catalogs from the low and 
high resolution runs overlap, the respective mass func¬ 
tions converge at ~ 5% level, while for the WDM-b case 
a scatter of ~ 20 % is found only in certain mass bins 
due to Poisson noise. Hence, this approach allows us to 
recover a numerically convergent mass function over a 
wider range of low masses compared to the aggressive 
Mim cut US] that discards all halos - genuine or other¬ 
wise - with mass M < My lm while still leaving the halo 
abundance uncorrected between M Vun and Mh m . This 
suggests that the halo selection criterion proposed here 
can efficiently remove spurious halos (i) independently of 
the resolution of the simulations, and (ii) independently 
of the specificities of the cosmological model. 

In Fig. [ 8 ] we plot the normalized probability density 
function of the spin parameter for all halos shown in 
Fig. a before and after removing unvirialized halos for 
two different cuts, 77 > 1.5 and 2 respectively. We see 
that removing unrelaxed halos from the numerical cat¬ 
alogs of the different models recovers spin distributions 
which are well approximated by a log-normal and are 
consistent with each other over the same range of values. 

In Fig. [9] we plot the normalized probability density 
function of the halo shape parameters after the cut 77 > 
1.5 for the same mass bins as shown in Fig. [5j Again we 
see that having removed unrelaxed halos gives consistent 
distributions of the shape parameters at low masses with 
minimal contribution from objects with extreme values 
of ellipticity and prolateness. 


V. CONCLUSIONS 

Unphysical groups of particles arise from the artificial 
fragmentation of the matter density field in N-body sim- 








FIG. 8: Probability density function of the spin parameter 
at z = 0 for ACDM (top left), LFDM (top right), WDM-a 
(bottom left) and WDM-b (bottom right) for all halos in the 
mass bins shown in Fig. [I] Green triangles correspond to the 
distributions inferred before removing unrelaxed halos, while 
the red circles and blue squares correspond to relaxed halos 
with 0 < r) < 1.5 and 0 < r/ < 2 respectively. The red and 
blue solid curves are the corresponding best-fitting log-normal 
functions. 

ulations of non-standard DM models characterized by a 
small scale cut-off in the linear matter power spectrum. 
These spurious halos contaminate numerical halo cata¬ 
logs altering the low-mass end of the halo mass function 
in the range of masses which are most relevant for ob¬ 
servational tests of DM scenarios. In the case of Warm 
Dark Matter models these artifacts manifest as an up¬ 
turn of the mass function at low masses. However, the 
absence of an upturn does not imply that the DM models 
with a softer suppression of power at small scales are ex¬ 
empt from spurious halo contamination, rather they are 
subdominant. 

Here, we have argued that the study of the structural 
properties of halos is a better proxy to differentiate be¬ 
tween artificial group of particles and genuine halos. Us¬ 
ing a suite of high-resolution simulations we have shown 
that spurious halos have systematically larger values of 
the spin parameter, are less spherical, more elliptical or 
prolate and strongly deviate from virial equilibrium com¬ 
pared to genuine halos. All these features are highly cor¬ 
related for spurious halos and we find them even in the 
case of the LFDM scenario whose mass function does not 
exhibit any upturn at low masses. 

The strong correlation between the spin parameter and 
the deviation from virial equilibrium suggests that select¬ 


ing approximately well relaxed halos can remove the bulk 
of spurious halo contamination. We find that a mildly 



- 0.1 0 0.1 0 . 20.1 0 0.1 0 . 20.1 0 0.1 0 . 20.1 0 0.1 0.2 


FIG. 9: Probability density function of halo sphericity, ellip- 
ticity and prolateness (panels top to bottom) at z = 0 for 
ACDM, LFDM, WDM-a and WDM-b (panels left to right) 
respectively for halos with virial state parameter 0 < r/ < 1.5. 
We consider the same mass bins as in Fig. [5] 


conservative cut in the virial state parameter allows to 
recover corrected halo mass function. Most importantly, 
we show that this approach gives convergent results when 
applied to halo catalogs computed from N-body simula¬ 
tions with different mass resolution. 

The advantage of this approach is twofold. One, it 
allows to estimate the correct halo abundance over a 
wider range of masses than a simple mass cut prescrip¬ 
tion which depends on the resolution of the simulation 
and the characteristics of the simulated cosmology. Two, 
as the halo selection is based on a physical criterion it 
is also applicable to cosmological models where the halo 
mass function lacks a visible upturn due to the spurious 
halos being subdominant at low masses. 
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